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Abstract —In this paper, the variable wind power is incorpo¬ 
rated into the dynamic model for long-term stability analysis. A 
theory-based method is proposed for power systems with wind 
power to conduct long-term stability analysis, which is able to 
provide accurate stability assessments with fast simulation speed. 
Particularly, the theoretical foundation for the proposed approx¬ 
imation approach is presented. The accuracy and efficiency of 
the method are illustrated by several numerical examples. 

Index Terms —wind power, stochastic differential equations, 
sufficient conditions, power system dynamics, power system 
stability 

I. Introduction 

Steady increase in load demands and aging transmission 
networks have pushed many power systems closer to their 
stability limit 12. Previous practical experience shows that 
the system may remain stable in the short-term time scale 
(say 0-30s after the contingency) while becomes unstable in 
the long-term time scale (say 30s-several minutes after the 
contingency) due to insufficient reactive power support or 
poor control schemes 12- 0- The focus of long-term stability 
is on slower and longer duration phenomena after generator 
dynamics damp out |4| - l6l . In this time scale, two-time scale 
decomposition of the power system dynamic model can be 
employed f2) j3). Particularly, the dynamics can be classified 
into short-term dynamics and long-term dynamics based on the 
time frames of function after contingencies. As the penetration 
of wind power continues to grow, the long-term stability of a 
system may be further affected by the volatile nature of wind 
power and the distinct dynamic response characteristics of 
wind turbines HI. Long-term stability analysis deserves more 
attention to better ensure the secure operation of power grids. 
It is also crucial to study the impacts of wind power on long¬ 
term stability. 

In previous literature, impacts of wind power on long-term 
stability have been discussed in 0 Col, which focus on 
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the effects of different control schemes of wind generators. 
Stability impacts of wind power on transient stability and small 
signal stability have been discussed in EQ- 113. In all above 
studies, power systems are formulated in ordinary differential 
equations (ODEs) or differential algebraic equations (DAEs), 
and wind speed is assumed to be constant without considering 
the statistic properties of wind. To study the effects of the 
variability of wind power, stochastic differential equations 
(SDEs) are employed in 133- ED for transient stability and 
small-signal stability. In those studies, the variations induced 
by wind speed are simply modelled as white noise perturba¬ 
tions on power injections without characterizing the statistic 
properties of wind speed. In ED, a set-theoretic method 
is proposed to assess the effect of variability of renewable 
energies on short-term dynamics under small perturbations. 
Regarding long-term stability study of power systems with 
wind power, the challenges include the stochastic character¬ 
ization of wind power and the resulting high computational 
burden in stability assessments, which have not been well 
addressed. Specifically, introduction of the randomness will 
clearly increase the computational burden, making the stability 
assessments of stochastic systems time consuming. 

In this paper, the statistic properties of wind power are 
well characterized and incorporated in stability analysis; a 
method is proposed with its theoretical foundation for long¬ 
term stability analysis, which is able to reduce the compu¬ 
tational burden arising from the randomness. Particularly, a 
SDE-based model E2 is utilized to characterize the Weibull- 
distributed wind speed, which is further incorporated into the 
complete dynamic model described in SDEs. Under this SDE 
formulation, a theory-based method, which is to approximate 
the stochastic model by a deterministic model, is proposed to 
perform long-term stability analysis. The theoretical founda¬ 
tion for the method is also developed under which the accuracy 
of the deterministic model is guaranteed. Compared to the 
deterministic models using constant wind speeds, the proposed 
deterministic model reflects the variable nature of wind power 
by providing correct stability assessments for the stochastic 
model. Compared to the stochastic model characterizing the 
statistic properties of wind speed, the proposed deterministic 
model takes much less time in time domain simulation. 

The rest of the paper is organized as follows. Section [II] 
briefly reviews power system models. Section [III] introduces 
the preliminaries about the wind speed model by SDEs and 
singular perturbation method. Section IIVI presents the for- 
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mulation of the long-term stability model with wind power 
based on SDEs. Afterwards, an analytical method with its 
theoretical foundation is proposed in Section [V] for long-term 
stability analysis of power systems with wind power. Several 
numerical examples are given in Section [Vi] to illustrate the 
feasibility and efficiency of the proposed method. Conclusions 
and perspectives are stated in Section IvTH 

II. Power System Models 

The deterministic power system long-term stability model, 
i.e., complete dynamic model, for simulating system dynamic 
response relative to a disturbance can be described as: 


4 

= eh c (z c ,x,y,z d ) 

( 1 ) 

X 

= f(zc,x,y,z d ) 

( 2 ) 

0 

= g{z c , x, y, z d ) 

(3) 

Zd(k) 

= h d {z c ,x,y,z d {k- 1 )) 

(4) 


Equation © describes long-term dynamics including ex¬ 
ponential recovery loads, turbine governors (TGs) and over 
excitation limiters (OXLs), and Eqn © describes the inter¬ 
nal dynamics of devices such as generators, their automatic 
voltage regulators (AVRs), certain loads such as induction 
motors, and other dynamically modeled components. Equation 
(0 describes the electrical transmission system and the internal 
static behaviors of passive devices, and Eqn © describes long¬ 
term discrete events like load tap changers (LTCs) and shunt 
switchings. h c , / and g are continuous functions, and vectors 

z c £ x £ R nx and y £ R ra « are the corresponding 

long-term/slow state variables, short-term/fast state variables, 
and algebraic variables, z d £ 1"^ are termed as long-term 
discrete variables whose transitions from z d (k — 1) to Zdik ) 
depend on system trajectories and occur at distinct times tk 
where k = 1,2,3Besides, 1/e can be regarded as the 
maximum time constant among devices. 

Load models play an important role in long-term stability 
analysis, and generally load dynamics appear in Eqn © as: 

p = x p /T p + p t (5) 

q = x q /Tq + q t (6) 

where x p and x q are internal state variables associated with 
generic load dynamics. They are described in Eqn © as: 

Xp — Xp/Tp T p s pi (7) 

Xq — X q /T q “t“ Qs qt (8) 

where p s and p t are the static and transient real power 
absorptions; similar definition for q s and q t . p° L and (// are 
PQ load power from power flow solutions a Ea. 

The tap-changing logic of LTC at time instant tk shows up 
in Eqn ©, and is given as follows: 

{ TO*, + Am if v > vq + d and mk < m max 
mk — Am if v < vq — d and m*, > rn mm 

mk otherwise 

(9) 

where v is the controlled voltage of LTC, vo is the reference 
voltage, d is half the LTC dead-hand, and rn max and m mm 
are the upper and lower tap limits Ei a Ei. 


III. Mathematical Preliminaries 

In this section, the method to model the Weibull-distributed 
wind speed by SDE proposed in ff2Tl j22l is briefly reviewed. 
Singular perturbation theory for ODE ED is also reviewed. 


A. Modelling Wind Speed by SDE 

Two continuous wind speed models based on SDE have 
been developed in ED ED- The developed models include 
but not limited to the autocorrelated Weibull distributed wind 
speed models. In this paper, model I in ED is applied to 
formulate wind speed. Briefly speaking, memoryless transfor¬ 
mation of Ornstein-Uhlenbeck process is used to obtain the 
Weibull distributed process. 

First, consider the following stochastic differential equation 
with given initial condition p(0) ~ A/”(0,/3 2 /2a): 


f] = —ap + P £ 


GO) 


where £(f) = dn d ^ is white noise process and W(t) rep¬ 
resents a standard Wiener process. The stochastic process 
r/(t) is an Ornstein-Uhlenbeck/Gauss-Markov process which 
is stationary, Gaussian and Markovian with the statistic prop¬ 
erties E [r)(t)} = 0, Var [rj(t)] = /3 2 /2a, and autocorrelation 
Aut = 

Then we apply a memoryless transformation to rj(t) as 
follows: 


l(t) = F w \<f>( 


v(.t) 

P/\[W) 


)) 


( 11 ) 


where F w is the Weibull cumulative distribution function: 

F w (u) = l-e (u/A) \ Vit > 0 ( 12 ) 


with the shape parameter k > 0 and the scale parameter A > 0, 
and $ is the Gaussian cumulative distribution function: 


,,u — Ektl. I., .. u — Efotl .. 

= 2 1 + erf ( /ov r 1 )} ’ 
y/Var|wJ z y^VarpiJ 


Vit £ 


(13) 


The resulting process 7 (t) is a Weibull distributed stochastic 
process with the statistic properties E[ 7 (t)] = AT(1 + i) = 
p w , Var[ 7 (f)] = A 2 T(1 + |) - pi, and Aut[ 7 ^), 'y(t j )] « 
e -oc\tj-u\' ]q ote t jj at E[ 7 (t)], Var[ 7 (f)] and Aut[ 7 (tj), 7 (tj)] 
do not depend on the parameter /3, i.e., /3 in Eqn (flTTb does not 
affect the statistic properties of stochastic process 7 (t) lf2TI . 
Since wind speed is autocorrelated in the range of hours, a 
depending on the decay rate of autocorrelation is usually small 
in the time frame of stability analysis. 

Note that although we consider the Weibull-distributed wind 
speed model in this paper, Eqn (ITOll -dTTIi can be used to 
formulate other probability distributions by applying different 
F w in Eqn dTTb . Refer to ED for more details. 


B. Singular Perturbation Theory for ODE 

We next consider the following general singular perturbed 
model (slow-fast model): 

£ e :i = f(z,x) z£K " 2 (14) 

ex = g(z,x) x £ l ” 1 
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where e is a small positive parameter, z is a vector of slow 
variables while x is a vector of fast variables. 

The slow model is obtained by setting e = 0 in (114b : 

E 0 : i = f{z,x ) z£l" ! (15) 

0 = g(z,x) i£R"‘ 

The algebraic equation 0 = g(z,x) constrains the slow 
dynamics to the constraint manifold defined by: 

T := {(z, x) G R"* x R n * : g(z, x) = 0} (16) 

The trajectory of model (Ti~5} starting at (zq. xq) is denoted by 
fo(f zo, xo), and the stability region is defined as: 

A 0 (z s ,:c s ) := {(z,x) G T : f 0 (t,z,x) -»• (z s ,x s ) as t -»• oo} 

The singular points of system <us> or singularity S is 
defined as: 

S := {(z,x) G T : det(d x g)(z,x) = 0} (17) 

where d x g is the Jacobian matrix of g with respect to x. 

Singular points can drastically influence the trajectories of 
the DAE system (Ti~5} . Typically, the singular set S' is a set of 
maximal dimension (n z — 1 ) embedded in T, and T is separated 
by S into open regions ||28l |29| . 

We next define the fast model, i.e., boundary layer model, 
associated with the singularly perturbed model. Define the fast 
time scale n = t/e. In this time scale, model ( 1 1 41 takes the 
form: 

dz 

n e : — = ef(z,x) ZGI"‘ (18) 

dK 

(IT 

-f = g(z, x) x G R™ x 

an 

Let e = 0, we obtain the fast model as follows: 

ID : ^ = 0 zGR" ! (19) 

UT 

= g{z, x) igR"- 

dn ; 

where z is frozen and treated as a parameter. The constraint 
manifold T is a set of equilibriums of models dT9l> . For each 
fixed z, a fast dynamical model (fl9t is defined. 

Definition 1: Stable Component of the Constraint Manifold / 
Uniformly Asymptotically Stable 
Assuming (z, x) ^ S, and x = j(z ) is an isolated root of 
equation: 

0 = g(z,x) ( 20 ) 

then x = j(z) is an equilibrium point of system (fl~9l) and is 
also termed as the slow manifold of system ( fl4l >. If x = j(z ) 
is an asymptotically stable equilibrium point of system ( IT9l ) 
for all zg 7, then we say that the equilibrium point j(z) is 
uniformly asymptotically stable , and the slow manifold x = 
j(z) is a stable component of the constraint manifold T. 

Under the following conditions, the slow model d~i~5} pro¬ 
vides good approximations of the singular perturbed model 
m Assume that the equilibrium point x = j(z) of the 
fast model is uniformly asymptotically stable in z G Z 
and the initial condition (zo,a;o) is in the stability region 
A t (zo-j(z 0 )) of the initial fast model with zo G Z, then 


solution of singular perturbed model (IT4} will approach j(z) 
in a time of order e|loge|, and the solution of model dT4t 
remains in an e-neighbourhood of the slow manifold j{z) 
[|23l li25l . Furthermore, there exists an invariant manifold 
x = j*(z, e) = j(z) + O(e) for sufficiently small e in z £ Z 
Erl lf30l . and the dynamics on invariant manifold is given by: 

z = f{j*(z,e),z) ( 21 ) 

which can be approximated by regular perturbation theory. 
Specifically, Eqn (fTil is reduced to slow model (ITtI i at e = 0. 

IV. SDE FORMULATION OF POWER SYSTEM MODELS 

When the long-term stability model is represented in r time 
scale, where r = te, and / refers to 4-, the long-term stability 
model of power system can be represented as: 


z 'c 

= h c (z c , x, y, Zd) 

( 22 ) 

ex' 

= f{z c ,x,y,z d ) 

(23) 

0 

= g(z c ,x,y,z d ) 

(24) 

z d (k) 

= h d (z c ,x,y, z d (k — 1 )) 

(25) 


It is possible to decouple the long-term stability model 
into two decoupled systems (I22t - (l24l > and ({25}. When discrete 
variables Zd jump from Zd(k— 1) to Zd(k), Zd are updated first 
according to d25t and then system d22} - (l24} works with fixed 
parameters Zd- In other words, the long-term stability model 
can be regarded as a collection of continuous DAE systems 
with fixed discrete parameters. In the remainder of this paper, 
we focus only on the continuous DAE system (l22} - (l24t . 

When incorporating wind power in the power system mod¬ 
els, we use Eqn (ITOt-dTIb to formulate wind speed. Assuming 
that there are n w wind sources in the system, and the wind 
model is described in r time scale, we have 

= - A Vw + <j£,= fwiVw) + (26) 

Vw = F- l {^{—j^=))=g w { rj w ) 

where rj w ,y w G R”“, A = diag[ai,..., q„J G 
in which oti is determined by autocorrelation 
property of each wind source, o is a small positive 
parameter, f Q £(s)ds is a n^-dimensional Wiener pro¬ 
cess, F w — \F W {j]wi ) ) F w (T]w 2 ) 5 **-Fw (jlw nw )] 5 ^ — 

[®(ri Wl ),${r) W2 ), ...$(r) Wnw )] T , and g w : R"“ i-aR”"". 

Note that a does not affect statistic properties of y w , 
which has the same reason that (3 does not affect the statistic 
properties of 7 (t) as shown in the last paragraph of Section 
IIII-AI This property of the formulation will play an important 
role in the analysis of Section IY1 

With (|22} - (|2?1 > and ({26}, the stochastic long-term stability 
model can be represented as: 

z' = h c (z c , x 1 y, z d ) (27) 

ex' = f(z c , x, y, z d ) + crB^ 

0 = g{z c ,x,y,z d ) 

where x = [x T ,g'^,] T G R n *+ n ” = R n ®, y = [y T ,y'^] T G 
R n » +n » = R n s, B = 
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Under the following assumption, the stochastic DAE system 
( OH can be simplified as a stochastic ODE system. 

Assumption 1. The trajectory of the deterministic long-term 
stability model OH doesn’t meet singularity. 

z' c = h c (z c ,x,y,z d ) (28) 

ex' = f(z c ,x,y,Zd) 

0 = g{z c ,x,y,z d ) 

Under Assumption 1, y can be expressed by (z c , x, z d ) from 
the algebraic equation as y = l{z c ,x, z d )- Substituting the 
expression back to model (f27l) . we have: 

z'c = H c (z c ,x,z d ) (29) 

ex' = F(z c , x, Zd) + <tB£ 

or written in the form: 

dz c = H c (z c , x, Zd)dr (30) 

dx = -F(z c ,x, z d )dr + -^=BdW 

e V e 

Besides, system (OH can be written as: 

dz c = H c {z c , x, Zd)dT (31) 

dx = —F{z c , x, z d )dT 

e 

In the next section, we will provide a theoretical foundation 
for approximating (OH by ( OH in long-term stability analysis. 
The deterministic model (OH can provide correct stability 
analysis results for the stochastic model (OH provided some 
mild conditions are satisfied. Compared to the stochastic model 
(Pol l, the deterministic model (OH takes less time in time 
domain simulation. Compared to other deterministic models 
using constant wind speed, the deterministic model (OH re¬ 
flects the variability of the wind speed by providing accurate 
stability assessments for (OH . 

V. THEORETICAL BASIS 

In this section, we apply the singular perturbation theory 
for SDEs proposed in j24l - |s26), to develop the theoretical 
foundation for approximating (Ool> by (OH . The main results 
of lf24i - lf26i are briefly reviewed in Appendix [A] Note that 
we have formulated the stochastic long-term stability model 
(Pol l into the conventional form (OH considered in singular 
perturbation theory for SDE. Specifically, Q = 0. We further 
assume the counter-part of Assumption 2.1 (see Appendix [A} 
in this specific formulation is satisfied. Note that differentia¬ 
bility and uniform boundedness are usually satisfied in real 
physical systems, and the assumption is reasonable. 

We denote the slow manifold of system (ED as x = 
h{z c ,Zd)- If x = li(z c , Zd) is a stable component of the 
constraint manifold when z c £ D~ c C R'~' : , then there exists 
an invariant manifold of system (OH that x = l\{z c ,Zd,e) = 
li (z c , z d ) + O(e). We further define the ellipsoidal layer N ( h ) 
surrounding the invariant manifold as below: 

N(h ) = {{z c ,x, z d ) : {{x - lf(z c ,z d ,e)), 

Li(z c , z d , e)~*\x - l*(z c , z d , e))) < h 2 }( 32) 


where () denotes the operator of inner product. Note that 
the matrix L\{z c ,z d ,e) describing the cross section of N[h) 
is well defined as shown in Appendix [B] In addition, an 
illustration of N(h) is shown in Fig. Q] 


Xi 



Fig. 1. Illustration of N(h) in which case n z = 1, n T ~2. N(h) is a 
ellipsoidal layer surrounding the invariant manifold l*(z c , 

Theorem 1 shows that under some mild conditions, the 
sample paths (i.e., trajectories) of the stochastic long-term 
stability model (QOt are concentrated in the ellipsoidal layer 
N(h) surrounding the invariant manifold l\{z c , z d ,e) of the 
deterministic long-term stability model ( OH . 

Theorem 1: Consider the systems (POt-tPH. for some fixed 
eo > 0, ho > 0, there exist 5o > 0, a time tt of order e|log/i| 
such that if the following two conditions are satisfied: 

i. The slow manifold x = h(z c ,Zd) where z c £ 
D Zc C R n * c is a stable component of the constraint 
manifold; 

ii. The initial condition (z c (0), x(0),Zd) of the stochas¬ 
tic long-term stability model (QOt where 2 C (0) £ D Zc 
satisfies (xr c (0), x(0), z d ) £ N(6q), 

then for r > tt, the sample path (z c (t), x(t), z d ) of the 
stochastic long-term stability model (QOt satisfies the following 
relation 

P{3s £ [tt,t] : ( z c ( r ), x ( r ), Zd ) $■ N(h)} 

< C nzc!nx (r,e)e^- (1 - OW -°^ ) (33) 

for all e < eo, h < ho, where the coefficient C nzci n x (T,e ) = 
[C n *c + h~ nx ]( 1 + is linear time-dependant. 

Proof: It can be directly deducted from Theorem 3 stated in 
Appendix lAl hence is omitted. 

Note that the detailed expression of C Uz jHx (r, e) is omitted 
here for brevity, which can be obtained directly from Theorem 
2.4 in |[25l . One key point is that the coefficient is independent 
of h and a, so the probability of the sample path to leave the 
layer N(h) decays exponentially as h increases. 

One important implication of the theorem is that as soon 
as we take h slightly larger than er, say h = O(ojlogaj), the 
right hand side of Eqn (OH becomes very small unless we 
wait for long time spans. That means the sample paths of the 
stochastic long-term stability model (QOt are concentrated in 
N(a). Furthermore, if h cr in Eqn (OH . sample pathes of 
the stochastic model unlikely leave N[h) as long as slow dy¬ 
namics permit. Hence, we can explore the trajectory relations 
between the stochastic model (OH and the deterministic model 
(OH without concerning probabilities. On the other hand, the 
condition that h cr can be easily satisfied in this SDE 
formulation of the power system model. As stated before, <j 
doesn’t affect the statistic properties of wind speed, thus we 
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can choose a as small as desired such that any satisfactory 
depth h of the layer satisfies er. 

Next we study the trajectory relations between the stochastic 
long-term stability model (l30l > and deterministic long-term 
stability model (ED under the condition that h cr, i.e., 
the depth h of layer N(h) is much larger than the small 
positive parameter a associated with wind speeds. According 
to singular perturbation theorems E3 ED and the theoretical 
foundation for the Quasi Steady-State (QSS) model proposed 
in EH, we have that if the slow manifold h(z c ,Zd) of the 
deterministic long-term stability model is a stable component 
of the constraint manifold, and the initial point on trajectory 
of the deterministic long-term stability model lies inside the 
stability region of the initial short-term stability model, then 
the trajectory of the long-term stability model will approach 
the slow manifold in a time of order n = e|loge|, and moves 
along the invariant manifold l*(z c , Zd, e). Refer to ||32| for 
more details. 

We denote the solution of stochastic long-term stabil¬ 
ity model gsd as (z c (T),x(T),Zd), and denote the so¬ 
lution of deterministic long-term stability model ED as 
(z c d{t),xd(t), Zd)- Then the following theorem describes the 
trajectory relations between the two models. 

Theorem 2 (Trajectory Relationship): Assuming that h cr, 
consider the systems ©-ED- for some fixed eo £ (0, h), 
there exist 5 0 > 0, a time t\ of order e\\ogh\, and r 2 > T\ 
such that if the following conditions are satisfied: 

i. The slow manifold x = li(z c ,Zd) where z c £ 
D Zc C l" Zc is a stable component of the constraint 
manifold; 

ii. The initial condition of stochastic system 

(2 c (0),i(0),Zd) where z c ( 0) £ D~ c satisfies 

(2 c (0),z(0),z d ) e N(S 0 ); 

iii. The initial condition (z c (0), xd( 0), Zd) of the de¬ 
terministic long-term stability model (OTT) is inside 
the stability region of the initial short-term stability 
model, 

then for r £ [ 71 , 72 ], the following relations hold: 

\x(t)-x d (t)\ = O(o) (34) 

\z c (t) - z c d(t)\ = (35) 

for all e £ (0, eo). 

Proof: See Appendix [C] 

Theorem 2 indicates that trajectory of the stochastic long¬ 
term stability model is concentrated in N(<r), and can be 
approximated by trajectory of the deterministic long-term 
stability model as shown in Fig. [2] provided that the required 
conditions are satisfied. On the other hand, these conditions are 
mild and are usually satisfied in practical applications, hence, 
the deterministic long-term stability model can be applied for 
stability analysis in most situations. 

VI. Numerical Illustration 

In this section, three examples are given to show that 
the deterministic long-term stability model can provide ac¬ 
curate approximations and correct stability assessment for 
the stochastic long-term stability model, provided that all 



Fig. 2. The trajectory z c (0), x(0), Zd) of the stochastic long-term 

stability model is concentrated in 7V(oj, and can be approximated by 
trajectory <f>i (r, z c (0), xd( 0), Zd) of the deterministic long-term stability 
model. 

conditions of Theorem 2 are satisfied. At the same time, the 
time consumed by the deterministic long-term stability model 
is much less than that consumed by the stochastic long-term 
stability model. All simulations were performed on PSAT-2.1.8 
|20| . Note that both the stochastic model and the deterministic 
model applied fixed time step size, whereas the efficiency 
of the deterministic model can be further improved by using 
larger or adaptive step size. 

A. Numerical Example 1 

The first example is a 14-bus system in which there are 
two wind sources and two doubly-fed induction generators 
(DFIGs) at Bus 2 and Bus 8 respectively. The wind speeds 
are characterized by Weibull distributed stochastic processes 
with parameters k = 1.51, A = 3.36 for wind source at Bus 
1 and k = 5, A = 20 for the other one. Both of them 
have autocorrelation parameter a = 0.2575/3600 ED- In 
addition, there are three synchronous generators (SYNs) which 
are controlled by AVRs, TGs and OXLs. And there is an 
exponential recovery load at Bus 5 and one continuous LTC 
between Bus 4-9 whose initial time delay is 20s. At Is, there 
is one line loss between Bus 10-9. 

One of the sample paths of the stochastic model took 
156.96s to simulate, while the deterministic model only took 
93.78s. The time required by one sample path has already 
been longer than that of the deterministic model, let alone the 
cases where several sample paths of the stochastic model are 
required to evaluate the stability of the system. The compar¬ 
ison between trajectories of the stochastic long-term stability 
model and those of the deterministic long-term stability model 
is shown in Fig. [3] Slow manifold of the deterministic model 
which is obtained from the Quasi Steady-State model El is 
also denoted. Trajectory of the stochastic long-term stability 
model is fluctuating around that of the deterministic long-term 
stability model due to variable wind speed, while it always 
keeps a small distance to trajectory the deterministic long¬ 
term stability model. The stability assessments of both models 
are the same that the system is long-term stable. In this case, 
the deterministic model doesn’t meet singularity as seen from 
Fig- 0 and the slow manifold of the deterministic model is 
a stable component of the constraint manifold; trajectory of 
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the deterministic long-term stability model moves along the 
invariant manifold l*(z c , Zd, e) which is a e-neighbourhood of 
the slow manifold, thus we can apply the results of Theorem 
2. As a result, the deterministic long-term stability model 
provides accurate approximations for the stochastic long¬ 
term stability model with correct stability assessment that the 
system is long-term stable. 


„ °- 708 
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Fig. 3. The trajectory comparisons of the stochastic long-term stability model 
and deterministic long-term stability model. Since all conditions of Theorem 
2 are satisfied, the deterministic long-term stability model provides accurate 
approximations in terms of trajectories and stability assessment. 


B. Numerical Example 11 

The second example is also a 14-bus system. There is one 
wind source and a corresponding DFIG at Bus 2. The Weibull 
distributed wind speed is characterized by parameters k = 
1.51, A = 3.36 and a = 0.2575/3600. In addition, there are 
four generators which are controlled by AVRs, TGs and OXLs. 
And there are three exponential recovery loads at Bus 9, Bus 
10 and Bus 14 respectively. There are two continuous LTCs 
between Bus 4-9 and Bus 14-13, whose initial time delays are 
30s. At Is, there are three line losses between Bus 6-13, Bus 
7-9 and Bus 6-11. 

In this case, both the stochastic long-term stability model 
and the deterministic long-term stability model have voltage 
collapse, and the deterministic long-term stability model suc¬ 
cessfully captures the voltage instability. One of the sample 
paths of the stochastic model took 120.86s to simulate, while 
the deterministic model took 90.64s. The comparison between 
the trajectories of the stochastic long-term stability model 
and those of the deterministic long-term stability model as 
well as the slow manifold are shown in Fig. [4] Similarly, as 
all conditions of Theorem 2 are satisfied, the deterministic 
long-term stability model provides accurate approximations for 
the stochastic long-term stability model with correct stability 
assessment. 

Physically speaking, as LTCs try to recover voltages at Bus 
9 and Bus 13, OXLs at Generator 2 and Generator 4 are 
activated. LTCs require more power support from generators, 
however, OXLs protect generators from overheating, thus 
restrict power output. As the system can not meet the reactive 
power demand imposed by LTCs, voltages finally collapse. 
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Fig. 4. The trajectory comparisons of the stochastic long-term stability 
model and the deterministic long-term stability model. Since all conditions of 
Theorem 2 are satisfied, the deterministic long-term stability model provides 
accurate approximations in terms of trajectories and stability assessment that 
the system is long-term unstable. 


C. Numerical Example III 

The last example is a larger system which has 145 buses 
ll33l . There are 6 DFIGs driven by different wind sources. 
Besides, there are 50 synchronous generators each of which 
is controlled by AVR, and Generator 10-Generator 20 are also 
controlled by TGs. There are OXLs at Generator 1-Generator 
6 whose initial time delays are 50s. Additionally, there are 
5 continuous LTCs between Bus 1-33, Bus 1-5, Bus 79-95, 
Bus 79-92, Bus 60-95 respectively. All LTCs have initial time 
delays of 50s. At 0.5s, there are three line losses including 
Bus 94-95, Bus 94-138 and Bus 95-138. 
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Fig. 5. The trajectory comparisons of the stochastic long-term stability model 
and the deterministic long-term stability model. The deterministic long-term 
stability model provides accurate approximations for both trajectories and 
stability assessment. 


The stochastic long-term stability model took 145.71s to 
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tic long-term stability model and those of the deterministic 
long-term stability model is shown in Fig. 0 Slow manifold 
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of the deterministic model is also labelled. The simulation 
results show that the deterministic long-term stability model 
is able to provide accurate approximations for the stochastic 
long-term stability model with correct stability assessment, as 
all conditions of Theorem 2 are satisfied. 

VII. Conclusions and Perspectives 

In this paper, we have incorporated the Weibull-distributed 
model of wind speed proposed in |2TI into the power system 
dynamic model for long-term stability analysis. Under the SDE 
formulation, we have theoretically shown that the stochastic 
model can be approximated by a corresponding determinis¬ 
tic model for stability analysis under mild conditions. The 
proposed deterministic model takes much less time than the 
stochastic model, while is able to provide accurate stability 
assessments for the stochastic model. 

This paper may represent the first attempt to apply the sin¬ 
gular perturbation method for stochastic differential equations 
in power system stability research. The proposed methodology 
seems promising for stability analysis as the analytical results 
developed are not limited by the statistic properties of wind 
power, specific system models, penetration level and so forth. 
More particularly, the proposed method can accommodate 
wind speed models with different distributions like Rayleigh 
distribution. Lognormal distribution. Gamma distribution, etc. 
As a systematic method to model power systems as SDEs 
has been proposed in El . we believe that the analytical 
results in this paper can be readily generalized, and thus 
provide a unified framework to conduct stability analysis for 
systems with uncertainties arising from load variations and 
other renewable energies. 

As the integration of wind power continues to grow, the 
concerns about its stability impacts outstand. Practically speak¬ 
ing, power system operations on one hand concern about the 
variability of wind power, whereas on the other hand regard 
the wind power as constant in stability analysis to reduce the 
computational cost. The theory-based method proposed in the 
paper takes the statistic properties of wind power into account, 
and meanwhile does not increase the computational burden. 
Specifically, the computational cost for the proposed method 
is the same as the model using constant wind speed, which 
has been widely applied |9)- 071 : the proposed deterministic 
model can also be easily incorporated in any dynamic model 
originally used in the time domain simulation. In addition, 
the analytical results in the paper provide an answer to the 
question that when the variability of wind power needs to 
be considered in long-term stability study, and when it does 
not. Those results will give a practical guidance for power 
system stability assessments if efficient numerical schemes 
can be developed to check the sufficient conditions presented 
in Theorem 2. For other future works, an extension of the 
proposed method to the stochastic long-term stability model 
with discrete variables is expected. In addition, the proposed 
method can be further improved by applying the QSS model 
0 0- As such, the QSS model is able to incorporate wind 
generation into the formulation, and thus provides a more 
efficient computational tool to conduct long-term stability 


analysis for power grids with wind power. This extension of 
the QSS model may be a valuable advance in dynamic stability 
analysis, considering the wide industry applications of the QSS 
model and the continuous increasing of wind penetration. The 
analytical results proposed in this paper will be the foundation 
for those exciting subsequent works. 

Appendix A 

SINGULAR PERTURBATION THEORY FOR SDE 

We consider the following stochastic system: 

dz = f(z,x)dt + dQ(z,x)dW z € HE (36) 
dx = -g(z,x)dt + —^=G(z,x)dW x £ R" x 

where W denotes a fc-dimensional standard Wiener process, 
and Q £ R” zXfe , G £ R n * xk are matrix-valued functions. We 
are interested in the case that d does not dominate a, i.e., 
d = po where p is uniformly bounded above in e. 

Assumption 2. 

2.1 In some open connected set D C R" z x M™*, / £ 
C 2 (D,R n *), g £ C 2 (D,W nx ) and Q £ C 1 {D, R nzXfe ), G £ 
C 1 (D,R nxXk ). f, g, Q, G and all their partial derivatives up 
to order 2, respectively 1, are uniformly bounded in norm in 

D. 

2.2 The slow manifold x = j(z) where z £ Dq C R nz of 
deterministic system (ITdl) is a stable component of constraint 
manifold T. 

2.3 The diffusion matrix G(z,j(z))G(z,j(z)) T is uni¬ 
formly positive definite in z £ Dq. 

Theorem 2.4 in l25l which is reviewed as Theorem 3 as be¬ 
low states that under Assumption 2, the sample path of system 
(l36l) are concentrated in an ellipsoidal layer surrounding the 
invariant manifold j*(z,e) of system (fl4l) : 

M{h) = {{z, x) : ((x-j*(z,e)),J*(z,e)~ 1 (x-j*(z,e))) < h 2 } 

2 (37) 

up to time t, with a probability roughly like l—(t 2 /e)e~ h / 2cr . 

In particular, J*(z,e) describing the cross section of the 
ellipsoid is well defined under Assumption 2. Refer to Section 

2.3 of l25l or Chapter 5.1.1 of ll24ll for more details. 

The following theorem shows that the sample paths of the 
stochastic system (f36l) are concentrated in a a-neighborhood 
of invariant manifold j*(z , e) of the deterministic system (IT4l) . 

Theorem 3 4251/ : Under Assumption 2, for some fixed 
eo > 0, ho > 0, there exist <5o > 0, a time t\ of order 
e|log/i| such that whenever S < So, if an initial condition of 
stochastic system (l36t (z(0),a:(0)) where c(0) £ Dq satisfies 
(z(0), x(0)) £ M(5), then for r > n, the following relation 

P{3s £ [tt,t] : {z(t),x(t)) M(h)} 

<Cn„»,(T,e)e# (1 -° ( ‘ ) - 0(£)) (38) 

holds for all e < eo, h < ho, where coefficient C nz , nx (T, e) = 
[C nz + h~ Ux ]( 1 + \) is linear time-dependant. 

A matching lower bound given in Theorem 2.4 of l l25l 
shows that the above bound is sharp in the sense that it 
captures the correct behaviour of the probability. 

Remark T. It is shown in Theorem 5.1.17 in lf24ft (Remark 
2.7 and Lemma 3.4 in ESI) that the stochastic singular 


perturbed system (l36l ) can be approximated by the reduced 
deterministic system: 

dz = f(z,j*(z,e))dt (39) 

and the deviation between the solution of ( l39l > and that of 
(l36l t is of order up to Lyapunov time of system ( l39l ). 
That means zd(C) provides 0(c r-y/e) approximation for z(t), 
i.e., \z(t) - ZD(t) | = 0(cry/c). 

Appendix B 

CROSS SECTION OF N(h) 

Following the deduction of .J*(z,e) in Section 2.3 of J25] 
or Chapter 5.1.1 of |[24l . we can readily obtain the following 
results. The cross section L\(z c , Zd, e) of N(h) is the solution 
of the following slow-fast system: 

= H c {z,l*{z c ,z d ,e)) (40) 

eLi = A(z c , Zd, e)L\ + L\A(z c , Zd, e) T + 

where 

A(z c , Zdi r) — d s F(z c , li(z c , Zd, e)) 

edzJi (Zci Zdi C)dxH c (z c: l^ ( z c , Zd, c)) 

and k > 0 is a sufficiently small positive parameter. 

Note that the exact expression of L\(z c , Zd,e) is not of 
interest and thus does not need to be calculated, and we 
only need to know that L\(z c , Zd,e) is well defined in this 
formulation. 


nln x 0 

0 In w 


(41) 


Appendix C 
PROOF OF THEOREM 2 

Proof: According to Theorem 2, if the conditions that h^> a 
and (i)-(ii) in Theorem 2 are satisfied, then after a time 
of order e|log/i|, sample paths of the stochastic long-term 
stability model (l30l > are concentrated in N(a) which is a a- 
neighborhood of the invariant manifold l\(z c ,Zd,Cj. Particu¬ 
larly, the slow dynamics of the stochastic long-term stability 
model m are in the ()(o\j A) neighborhood of the invariant 
manifold l\(z c ,Zd,e) as explained in Remark 1. Hence, Eqn 
(f34b and Eqn d35l > hold only if trajectory of the deterministic 
long-term stability model moves along the invariant manifold 
l* (z c , Zd. e). On the other hand, from singular perturbation 
theory, we have that if the slow manifold x = h(z c ,Zd) is a 
stable component of the constraint manifold, then the invariant 
manifold l*(z c , Zd, e) = li(z c , Zd)+0(e) exists for sufficiently 
small e. Moreover, if one more condition that the initial point 
of the deterministic long-term stability model lies inside the 
stability region of the initial short-term stability model is 
also satisfied, then there exists eo such that trajectory of the 
deterministic long-term stability model (PH moves along the 
invariant manifold l*(z c ,Zd,e) for all e £ (0, eo), r £ [ti,t 2 ] 
|[23l |32 ). Thus Eqn P4l i and Eqn P5l) hold for all e £ (0, eo), 
t £ [tt , t 2 ] . This completes the proof. 
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